Hola, soy Aitor Gonzalez, Yo aprendo Tidymodels para que te sea más fácil de aprender.
Soy data scientist o lo que es lo mismo, estadístico que hace más cosas aparte de una regresión logística. Actualmente trabajo en la unidad de facultad de medicina y biológia de la universidad autónoma de Barcelona.
Este taller contiene mucho texto. Y suele requerir de un par de veces mirárselo para pillarle todo el tranquillo. No pasa nada si no te enteras a la primera. “Dentro de poco” estará subido a Youtube en mi canal Reestimando.
La historia comienza con este señor de la foto de la ziquierda. Él es Max Khun. Parece un señor muy majo o muy travieso. Yo personalmemte creo más lo segundo viendo su creación.
Comienza por aquí por que este señor fue el desarrollador de la Librería Caret.Caret fue una librería que comenzó en 2008 y que asentó las bases del modelado estadístico en R. Desde su creación, el paquete ha sido citado en más de 7000 publicaciones académicas, teniendo un enorme éxito en la comunidad.
El paquete fue un éxito debido a que antes de él, no había una forma de unificar sintaxis entre los principales paquetes de modelado. Cada modelo estaba en un paquete y todos diferentes entre si, un caos para organizar un flujo de trabajo. Caret permitía utilizar con una sintaxis común los modelos de varias librerías. Además también traía varias funciones para hacer automáticamente cosas que son necesarias en todo proceso de modelado, como una función que haga una partición de datos en un conjunto de entrenamiento y uno de prueba (training y test).
ELIPSIS en 2012 aparece la el entorno tidyverse de la mano de Hadley Whickam, la otra súper estrella de R en los últimos años. Tidyverse proporciona un montón de herramientas para poder lidiar con facilidad con la manipulación de datos, conviertiéndose en un packate súper popular no solo ne R, sino en toda la ciencia datos. El principal punto es la facilidad de sintaxis para todas sus funciones, haciendo que los usuarios tuvieran una lectura de código muy fácil y fluida, enfocandose 100% en la necesidad y no en la programación.
Encantado con un framework en el que toda la sintaxis está ordenada, en 2019, Max se une al tidyverso, y comienza a crear Tidymodels cooperando con el equipo de Whickam. El objetivo era convertir Careta algo más modularizado y que permitiera centrar en el modelado en si mismo, más que infinitas líneas de sintaxis para poder crear modelos de machine learning.
A fecha de este taller, noviembre de 2023, hay otra desarrolladora principal de Tidymodels además de MAx. Julia Silge.
Esta mujer, no solo es un as de la programación, sino que además se dedica a expandir la filosofía Tidy a través de sus redes. Su contenido vale oro y diría que es imperdible si quieres no solo aprender más de Tidymodels, sino de las posibilidades del aprendizaje automático en general.
- https://www.youtube.com/@JuliaSilge
- https://juliasilge.com/blog/
Entre algunos de sus trabajos más llamativos están uno que me encanta que es Topic modeling for #TidyTuesday Taylor Swift lyrics es decir, “modelizando los temas que tratan las canciones de taylor swift” https://juliasilge.com/blog/taylor-swift/
Imaginemos que queremos hacer una regresión lineal. Cogeremos unos datos cualquiera, lo de mtcars.supongamos que queremos hacer una regresión logítica, la sintaxis para eso sería:
Como se pude ver en el ejemplo, la regresión fuera de tidymodels la hacemos con la librería glm. En Tidymodels esto se especifica con la función set_engine.
Ahora, veamos que pasaría si queremos hacerla con las liberría h2o
Se va pillando la idea, ¿no? Tidymodels permite que con cambios mínimos de código puedas cambiar completamente tu flujo de trabajo. Si queremos cambiar un “motor” de modelo, solo tenemos que especificarlo en un solo cambio, y no tener que reescribir todo el código de nuevo. Esto es clave ya que permite jugar mucho a la hora de probar muchos tipos de modelos y ampliar las miras sin aumentar los problemas de código.
Evidentemente, esto no es ni la superficie de Tidymodels, aquí solo estamos viendo la isla a lo lejos.
Somos unos expertos en aprendizaje automático y nos piden ayuda de un hospital 🏥
Al parecer, nuestra variable, Var_respuesta, es una enfermedad complicada de diagnosticar. Los médicos no son capaces de determinar si un paciente tendrá o no esa enfermedad al cabo de un año. Entonces han hecho un estudio caso-control en el que se miran unas cuantas variables; y al cabo de un año se da un diagnostico definitivo de var_respuesta. En este caso pues, nos interesa una modelización de tipo clasificatoria binaria, es decir de 2 niveles, o tiene la enfermedad o no la tiene. 👍/👎
Además de de nuestra variable Var_respuesta, tenemos otras 25 variables de diversos tipos, de cognición, bioquímicas, genéticas, etc. El dataframe es una simulación de uno real, e incluye valores perdidos y otras características a tener en cuenta en la fase de pre-procesamiento.
Insisto, el propósito es modelizar un clasifiación binaria
Antes de empezar podemos ver un poco su estructura.
# |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
# Librerias | Carga de datos ----
# |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
if(!require('pacman')){install.packages('pacman')}
pacman::p_load(
readxl, # Lectura de los datos
tidyverse , # Acceso al entorno de procesamiento "tidyverse"
tidymodels, # Acceso al entorno de modelado "tidymodels"
agua, # Modelaje de entorno h2o
themis, # Soporte de la librería "recipes" para añadir sobremuestreo
skimr, # Descriptiva rápida de datos
naniar,
knitr, # pPesta en html para este formato
)
datos <- read_xlsx('Data/datos.xlsx') %>%
filter(!is.na(Var_respuesta) ) %>%
mutate(Var_respuesta= as_factor(Var_respuesta))
datos %>%
head(.,10) %>%
knitr::kable()| Var_respuesta | Edad | Cognicion_01 | Cognicion_02 | Cognicion_03 | Cognicion_04 | Cognicion_05 | Bioquimica_01 | Bioquimica_02 | Bioquimica_03 | Bioquimica_04 | Bioquimica_05 | Bioquimica_06 | Bioquimica_07 | Bioquimica_08 | Bioquimica_09 | Bioquimica_10 | Bioquimica_11 | Escala_01 | Escala_02 | Escala_03 | Genetica_01 | Genetica_02 | Genetica_03 | Genetica_04 | Genetica_05 | Genetica_06 | Genetica_07 |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 14 | 101.82 | NA | NA | NA | NA | NA | NA | 5.65 | 13.65 | 125.71 | NA | NA | 84.29 | 1.48 | 30.35 | 35.68 | 0.00000 | 28.397013 | 7.363183 | NA | 0.9665678 | NA | NA | 1.2921746 | 0.2417377 | NA |
| 0 | 17 | 80.41 | 85.53 | 55.68 | 124.52 | 159.70 | 62.37 | NA | 5.65 | 17.83 | 154.53 | NA | NA | 99.42 | 1.34 | 33.99 | 33.01 | 138.73119 | 16.769985 | 17.764831 | 0.4550571 | NA | NA | NA | NA | NA | 0.0044324 |
| 0 | 14 | 91.00 | 79.76 | 94.28 | 140.68 | 100.07 | 86.33 | 26.24 | 3.18 | 13.09 | 135.88 | 81.53 | 5.75 | 109.45 | 0.95 | 30.44 | 34.21 | NA | 21.382220 | 7.061147 | NA | NA | NA | NA | NA | NA | NA |
| 0 | 21 | 86.99 | 91.33 | 93.65 | 115.00 | 123.01 | 66.84 | 10.23 | 3.48 | 11.97 | 132.15 | 90.66 | 2.58 | 68.60 | 1.04 | 30.10 | 34.00 | 48.83394 | 7.785552 | 20.442332 | NA | NA | NA | NA | NA | NA | NA |
| 1 | 33 | 78.83 | 104.19 | 42.68 | 112.97 | 179.16 | 59.51 | NA | 4.05 | 13.99 | 118.33 | 114.46 | 3.03 | 62.87 | 1.14 | 22.08 | 33.50 | 33.37826 | 28.165315 | 28.817018 | NA | 1.0690024 | 0.0339304 | 0.7820054 | 0.3230982 | NA | 0.3998781 |
| 1 | 27 | 93.54 | 81.16 | 89.32 | 115.83 | 65.07 | 93.17 | 32.35 | 3.50 | 17.95 | 132.61 | 67.61 | 4.10 | 76.30 | 1.11 | 36.50 | 32.24 | 33.45837 | 7.376139 | 22.080631 | 1.3081069 | 0.5763763 | NA | 0.5438819 | NA | 0.4705209 | NA |
| 0 | 16 | 114.25 | NA | NA | NA | NA | NA | 8.14 | 5.68 | 12.25 | 149.89 | NA | NA | 91.22 | 1.28 | 35.86 | 40.11 | NA | 18.489371 | 10.365972 | 1.5581016 | NA | NA | 1.1597700 | NA | 0.6242926 | 0.0579575 |
| 0 | 13 | 74.80 | NA | NA | NA | NA | NA | 93.32 | 4.77 | 8.88 | 164.33 | NA | 4.33 | 84.49 | 1.07 | 23.16 | 27.76 | NA | 18.749210 | 9.068513 | NA | NA | NA | NA | NA | NA | NA |
| 0 | 14 | 66.78 | NA | NA | NA | NA | NA | 54.33 | 2.56 | 16.07 | 161.64 | NA | NA | 99.47 | 0.92 | 30.14 | 34.27 | 133.14928 | 18.116769 | 31.680105 | 1.6247538 | NA | 0.9775707 | NA | 0.9699578 | NA | 1.1983120 |
| 1 | 17 | 84.13 | 100.23 | 69.62 | 137.84 | 91.24 | 78.76 | 28.80 | 4.22 | 14.38 | 151.67 | NA | 5.23 | 91.19 | 1.75 | 18.24 | 37.76 | 113.35402 | 20.327740 | 21.101902 | NA | NA | 1.1801013 | NA | NA | 1.0568568 | NA |
| Name | datos |
| Number of rows | 285 |
| Number of columns | 28 |
| _______________________ | |
| Column type frequency: | |
| factor | 1 |
| numeric | 27 |
| ________________________ | |
| Group variables | None |
Variable type: factor
| skim_variable | n_missing | complete_rate | ordered | n_unique | top_counts |
|---|---|---|---|---|---|
| Var_respuesta | 0 | 1 | FALSE | 2 | 0: 227, 1: 58 |
Variable type: numeric
| skim_variable | n_missing | complete_rate | mean | sd | p0 | p25 | p50 | p75 | p100 | hist |
|---|---|---|---|---|---|---|---|---|---|---|
| Edad | 9 | 0.97 | 22.33 | 6.47 | 8.00 | 17.00 | 21.00 | 27.00 | 39.00 | ▂▇▆▅▂ |
| Cognicion_01 | 0 | 1.00 | 86.51 | 29.37 | 0.00 | 78.50 | 90.56 | 103.13 | 143.84 | ▁▁▅▇▂ |
| Cognicion_02 | 82 | 0.71 | 90.22 | 12.98 | 55.23 | 80.74 | 90.43 | 98.98 | 129.15 | ▁▆▇▃▁ |
| Cognicion_03 | 57 | 0.80 | 77.62 | 17.80 | 35.22 | 65.91 | 76.77 | 88.23 | 130.09 | ▂▇▇▃▁ |
| Cognicion_04 | 64 | 0.78 | 134.85 | 51.16 | 15.31 | 97.17 | 137.91 | 172.32 | 245.38 | ▂▅▇▆▂ |
| Cognicion_05 | 70 | 0.75 | 126.23 | 44.38 | 53.89 | 94.13 | 120.96 | 155.32 | 249.33 | ▆▇▆▂▁ |
| Bioquimica_01 | 54 | 0.81 | 75.94 | 14.24 | 42.80 | 66.92 | 74.56 | 84.50 | 118.52 | ▂▇▇▃▁ |
| Bioquimica_02 | 43 | 0.85 | 36.07 | 38.35 | 0.12 | 11.43 | 24.47 | 42.21 | 250.61 | ▇▂▁▁▁ |
| Bioquimica_03 | 13 | 0.95 | 4.05 | 1.64 | 1.27 | 2.90 | 3.67 | 4.74 | 12.68 | ▇▇▂▁▁ |
| Bioquimica_04 | 12 | 0.96 | 14.30 | 1.96 | 8.88 | 12.95 | 14.34 | 15.54 | 20.65 | ▁▆▇▃▁ |
| Bioquimica_05 | 18 | 0.94 | 140.83 | 15.09 | 99.88 | 129.72 | 141.38 | 152.12 | 187.37 | ▁▆▇▃▁ |
| Bioquimica_06 | 54 | 0.81 | 86.32 | 13.75 | 55.96 | 77.21 | 86.04 | 95.10 | 142.74 | ▃▇▅▁▁ |
| Bioquimica_07 | 34 | 0.88 | 3.58 | 1.15 | 0.83 | 3.08 | 3.74 | 4.30 | 6.56 | ▂▂▇▃▁ |
| Bioquimica_08 | 16 | 0.94 | 84.37 | 12.35 | 43.50 | 76.65 | 83.96 | 91.39 | 126.11 | ▁▃▇▃▁ |
| Bioquimica_09 | 57 | 0.80 | 1.13 | 0.29 | 0.51 | 0.91 | 1.12 | 1.33 | 2.09 | ▃▇▆▂▁ |
| Bioquimica_10 | 12 | 0.96 | 30.09 | 3.88 | 18.24 | 27.48 | 30.34 | 32.51 | 40.62 | ▁▃▇▆▁ |
| Bioquimica_11 | 14 | 0.95 | 33.84 | 3.66 | 21.90 | 31.56 | 33.71 | 36.08 | 43.57 | ▁▃▇▆▂ |
| Escala_01 | 32 | 0.89 | 103.18 | 122.93 | 0.00 | 21.23 | 54.69 | 137.77 | 892.17 | ▇▂▁▁▁ |
| Escala_02 | 0 | 1.00 | 18.50 | 8.08 | 6.17 | 11.19 | 18.43 | 23.58 | 41.45 | ▇▇▇▃▁ |
| Escala_03 | 0 | 1.00 | 18.65 | 8.24 | 5.79 | 12.32 | 17.76 | 24.00 | 46.96 | ▇▇▅▂▁ |
| Genetica_01 | 160 | 0.44 | 0.80 | 0.60 | 0.00 | 0.29 | 0.64 | 1.23 | 2.99 | ▇▅▃▁▁ |
| Genetica_02 | 170 | 0.40 | 0.78 | 0.55 | 0.01 | 0.33 | 0.73 | 1.05 | 2.27 | ▇▇▅▂▂ |
| Genetica_03 | 178 | 0.38 | 0.76 | 0.56 | 0.00 | 0.32 | 0.67 | 1.09 | 2.22 | ▇▇▅▂▂ |
| Genetica_04 | 168 | 0.41 | 0.84 | 0.55 | 0.01 | 0.42 | 0.70 | 1.20 | 2.60 | ▇▇▅▂▁ |
| Genetica_05 | 145 | 0.49 | 20.19 | 214.06 | 0.00 | 0.84 | 1.82 | 3.17 | 2534.81 | ▇▁▁▁▁ |
| Genetica_06 | 178 | 0.38 | 0.84 | 0.53 | 0.03 | 0.37 | 0.77 | 1.21 | 2.36 | ▇▇▆▂▁ |
| Genetica_07 | 137 | 0.52 | 0.83 | 0.59 | 0.00 | 0.32 | 0.72 | 1.25 | 2.35 | ▇▆▅▃▂ |
El flujo de trabajo para producir un modelo en tidymodels es el siguiente:
La idea detrás consiste en una sola sintaxis que resuma todo este flujo de trabajo.
Como hemos visto en el ejemplo anterior, con cambiar una sola línea línea, podemos cambiar el tipo de modelo. Y con entender unos pocos puntos clave de todos los nexos, podemos llevar a cabo flujos de trabajo muy eficientes con muy poco coste de mantenimiento y poco cambio de sintaxis.
Para la creación de este “cosmos”, Tidymodels utiliza las siguientes librerías:
Para comenzar a crear un modelo de aprendizaje automático es necesario tener datos, fundamental. Estos datos, sueen ser dificiles de conseguir, requieren tiempo, esfuerzo y dinero. Su manejo debe ser optimizado con el fin de abaratar los costes de recopialrlos y mejorar los resultados. Además tener en cuenta que el propósito de un modelo de aprendizaje automático en general es que se pueda usar para poder predecir que va a pasar con unos datos que no tenga, o que tenga en un futuro.
Aquí entre el punto del remuestreo. Crearemos un partición de la muestra en unos datos para “entrenar” el modelo y la otra para “evaluarlo”, para ver si tal y como ha quedado calibrado el modelo es útil o no.
Por supuesto, la calidad de esta partición estará ligada a la calidad de los datos. Si los datos recopilados no son suficientemente representativos, el modelo no predecirá bien cuando traigamos casos que no haya visto, aún incluso si la evaluación del modelo ha sido buena.
Básicamente queremos esto
La libreria dentro de tidymodels que nos dará los conjuntos de datos para hacer esto es rsample.
Algunas de las funciones importantes dentro de este conjunto:
# 🔴 La función principal de Rsample, hace una partición de un conjunto de 🔴
initial_split(datos, prop = 0.6 )## <Training/Testing/Total>
## <171/114/285>
## # Bootstrap sampling
## # A tibble: 5 × 2
## splits id
## <list> <chr>
## 1 <split [285/105]> Bootstrap1
## 2 <split [285/108]> Bootstrap2
## 3 <split [285/112]> Bootstrap3
## 4 <split [285/107]> Bootstrap4
## 5 <split [285/102]> Bootstrap5
## # Rolling origin forecast resampling
## # A tibble: 6 × 2
## splits id
## <list> <chr>
## 1 <split [20/20]> Slice1
## 2 <split [20/20]> Slice2
## 3 <split [20/20]> Slice3
## 4 <split [20/20]> Slice4
## 5 <split [20/20]> Slice5
## 6 <split [20/20]> Slice6
ggpubr::ggarrange(
rolling_origin(datos, initial = 20, assess = 20, skip = 40, cumulative = FALSE) %>%
tidy() %>%
ggplot(aes(x = Resample, y = factor(Row), fill = Data)) +
geom_tile() +
scale_fill_manual( values = c('cyan','orange'))
,
rolling_origin(datos, initial = 50, assess = 20, skip = 70, cumulative = FALSE) %>%
tidy() %>%
ggplot(aes(x = Resample, y = factor(Row), fill = Data)) +
geom_tile() +
scale_fill_manual( values = c('cyan','orange'))
,
rolling_origin(datos, initial = 20, assess = 20, skip = 40, cumulative = T) %>%
tidy() %>%
ggplot(aes(x = Resample, y = factor(Row), fill = Data)) +
geom_tile() +
scale_fill_manual( values = c('cyan','orange'))
,
rolling_origin(datos, initial = 50, assess = 50, skip = 70, cumulative = T) %>%
tidy() %>%
ggplot(aes(x = Resample, y = factor(Row), fill = Data)) +
geom_tile() +
scale_fill_manual( values = c('cyan','orange'))
,
common.legend = T, legend = 'bottom'
)Una de las cosas que pasan habitual a la hora de crear una partición de datos es que, al ser aleatorias, uno de los conjuntos contenga más datos de un tiempo de otro tipo.
Idealmente querríamos que el conjunto de entrenamiento y que el conjunto de evaluación contuvieran la misma proporción de las categorías.
Para ello, en varias de las funciones de rsample se cuenta con el argumento strata que permite hacer esto de form estratificada para una variable o variables de interés, es este caso, para nuestros datos y nuestra variable conjunta:
set.seed(4147)
# Splits ----
## Sin estratificar
Split_datos_no_strat <- initial_split(datos, prop = 0.70 )
# Put 3/4 of the data into the training set
# 🔴 Initial_split con datos estratificados 🔴
## Estratificados
Split_datos_strat <- initial_split(datos, prop = 0.70, strata = Var_respuesta)ggpubr::ggarrange(
ggpubr::ggarrange(
training(Split_datos_no_strat) %>%
group_by(Var_respuesta) %>%
summarise(recuento = n()) %>%
ungroup() %>%
mutate(prop_var_respuesta = recuento / sum(recuento) ) %>%
ggplot(aes(x = Var_respuesta, y = recuento, fill = Var_respuesta)) +
geom_bar(stat = "identity") +
geom_text(aes(label = scales::percent(prop_var_respuesta)), position = position_stack(vjust = 0.5), size = 4) +
scale_fill_manual( values=c("azure4", "azure3"))
,
testing(Split_datos_no_strat) %>%
group_by(Var_respuesta) %>%
summarise(recuento = n()) %>%
mutate(prop_var_respuesta = recuento / sum(recuento) ) %>%
ungroup() %>%
ggplot(aes(x = Var_respuesta, y = recuento, fill = Var_respuesta)) +
geom_bar(stat = "identity") +
geom_text(aes(label = scales::percent(prop_var_respuesta)), position = position_stack(vjust = 0.5), size = 4) +
scale_fill_manual( values=c("azure4", "azure3")),
common.legend = T, legend = 'bottom'
),
ggpubr::ggarrange(
training(Split_datos_strat) %>%
group_by(Var_respuesta) %>%
summarise(recuento = n()) %>%
ungroup() %>%
mutate(prop_var_respuesta = recuento / sum(recuento) ) %>%
ggplot(aes(x = Var_respuesta, y = recuento, fill = Var_respuesta)) +
geom_bar(stat = "identity") +
geom_text(aes(label = scales::percent(prop_var_respuesta)), position = position_stack(vjust = 0.5), size = 4)+
scale_fill_manual( values=c("#E86EB7", "#108064"))
,
testing(Split_datos_strat) %>%
group_by(Var_respuesta) %>%
summarise(recuento = n()) %>%
mutate(prop_var_respuesta = recuento / sum(recuento) ) %>%
ungroup() %>%
ggplot(aes(x = Var_respuesta, y = recuento, fill = Var_respuesta)) +
geom_bar(stat = "identity") +
geom_text(aes(label = scales::percent(prop_var_respuesta)), position = position_stack(vjust = 0.5), size = 4)+
scale_fill_manual( values=c("#E86EB7", "#108064")),
common.legend = T,
legend = 'bottom'
),
labels = c('No estratificados', 'Estratificados'), label.x = 0.1
) %>%
plot()| Aclaración: los colores elegidos para el gráfico son de mi escena preferida de Spiderman, into the spiderverse. dejad siempre “huevos de Pascua” en los trabajos, hace que las cosas sean más emocionantes :P |
Para esto usaremos las funciones training() y test()
# 🔴 Deshacemos las particiones con las funciones training y testing 🔴
Training_datos <- training(Split_datos_strat)
Testing_datos <- testing(Split_datos_strat)
knitr::kable(head(Training_datos,10))| Var_respuesta | Edad | Cognicion_01 | Cognicion_02 | Cognicion_03 | Cognicion_04 | Cognicion_05 | Bioquimica_01 | Bioquimica_02 | Bioquimica_03 | Bioquimica_04 | Bioquimica_05 | Bioquimica_06 | Bioquimica_07 | Bioquimica_08 | Bioquimica_09 | Bioquimica_10 | Bioquimica_11 | Escala_01 | Escala_02 | Escala_03 | Genetica_01 | Genetica_02 | Genetica_03 | Genetica_04 | Genetica_05 | Genetica_06 | Genetica_07 |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 14 | 66.78 | NA | NA | NA | NA | NA | 54.33 | 2.56 | 16.07 | 161.64 | NA | NA | 99.47 | 0.92 | 30.14 | 34.27 | 133.149284 | 18.116769 | 31.680105 | 1.6247538 | NA | 0.9775707 | NA | 0.9699578 | NA | 1.1983120 |
| 0 | 23 | 108.18 | 84.49 | 82.39 | 183.52 | 166.31 | 101.24 | 3.66 | 3.18 | 20.29 | 148.58 | 100.23 | 3.53 | 81.72 | 1.11 | 35.75 | 34.03 | 68.943373 | 12.183308 | 15.670177 | NA | 1.797438 | NA | NA | 0.1373531 | 1.1481697 | NA |
| 0 | 26 | 0.00 | NA | NA | NA | NA | NA | 129.41 | 5.99 | 14.98 | 145.04 | 94.53 | 4.20 | 78.76 | 0.87 | 32.57 | 26.89 | 2.898493 | 11.236520 | 6.628174 | NA | NA | NA | 1.6622978 | 2.6878823 | 1.6836672 | 0.1001570 |
| 0 | 14 | 52.88 | NA | NA | NA | NA | NA | NA | 2.06 | 13.98 | 132.94 | 105.75 | NA | 82.49 | NA | 34.73 | 29.83 | 99.600374 | 25.026820 | 12.949102 | 0.1247104 | NA | NA | 0.6977428 | 1.2486785 | 0.4195379 | 0.3760827 |
| 0 | 27 | 109.25 | 88.04 | 103.31 | 143.93 | 144.97 | 100.08 | 8.09 | 4.49 | 19.09 | 187.37 | 84.47 | 3.87 | 91.97 | 1.26 | 30.74 | 34.95 | 0.000000 | 16.798451 | 7.607729 | 0.8493442 | NA | 0.8781136 | 1.9962258 | NA | NA | NA |
| 0 | 30 | 105.60 | 104.09 | 86.45 | 149.99 | 105.02 | 93.65 | 7.20 | 4.51 | 15.65 | 145.00 | 90.13 | 3.57 | 90.67 | 1.05 | 28.98 | 31.50 | 51.087736 | 6.731747 | 10.894387 | NA | 1.349131 | 0.3128530 | NA | 1.9916249 | 1.2268437 | NA |
| 0 | 12 | 91.73 | NA | NA | NA | NA | NA | 30.29 | 3.15 | 13.94 | 134.92 | 103.40 | 3.97 | 96.92 | 1.27 | 21.42 | 29.26 | 160.752695 | 17.851526 | 18.519326 | 0.8739840 | NA | 0.5427362 | NA | 1.2463796 | 0.8887141 | 0.7917460 |
| 0 | 17 | 91.97 | 95.35 | 80.65 | 172.32 | 122.78 | 68.50 | 40.55 | 2.77 | 12.49 | 145.84 | 74.00 | 4.01 | 86.12 | 1.16 | 27.01 | 33.90 | 15.319511 | 14.058768 | 17.484487 | 0.2352628 | NA | 0.3039570 | NA | NA | NA | 0.9261596 |
| 0 | 14 | 69.65 | NA | NA | NA | NA | NA | 55.17 | 4.87 | 11.98 | 165.35 | 91.39 | 5.73 | 52.67 | 1.21 | 29.25 | 32.46 | 76.930164 | 22.999587 | 14.681546 | NA | NA | NA | NA | NA | NA | NA |
| 0 | 34 | 111.50 | 58.99 | 101.17 | 181.40 | 119.87 | 106.90 | 9.96 | 3.19 | 14.47 | 159.95 | 102.97 | 3.98 | 85.87 | 0.91 | 30.43 | 29.18 | 7.373436 | 8.195322 | 6.718927 | 1.2862183 | NA | NA | NA | NA | NA | 0.2123626 |
| Var_respuesta | Edad | Cognicion_01 | Cognicion_02 | Cognicion_03 | Cognicion_04 | Cognicion_05 | Bioquimica_01 | Bioquimica_02 | Bioquimica_03 | Bioquimica_04 | Bioquimica_05 | Bioquimica_06 | Bioquimica_07 | Bioquimica_08 | Bioquimica_09 | Bioquimica_10 | Bioquimica_11 | Escala_01 | Escala_02 | Escala_03 | Genetica_01 | Genetica_02 | Genetica_03 | Genetica_04 | Genetica_05 | Genetica_06 | Genetica_07 |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 14 | 101.82 | NA | NA | NA | NA | NA | NA | 5.65 | 13.65 | 125.71 | NA | NA | 84.29 | 1.48 | 30.35 | 35.68 | 0.00000 | 28.397013 | 7.363183 | NA | 0.9665678 | NA | NA | 1.292175 | 0.2417377 | NA |
| 0 | 17 | 80.41 | 85.53 | 55.68 | 124.52 | 159.70 | 62.37 | NA | 5.65 | 17.83 | 154.53 | NA | NA | 99.42 | 1.34 | 33.99 | 33.01 | 138.73119 | 16.769985 | 17.764831 | 0.4550571 | NA | NA | NA | NA | NA | 0.0044324 |
| 0 | 14 | 91.00 | 79.76 | 94.28 | 140.68 | 100.07 | 86.33 | 26.24 | 3.18 | 13.09 | 135.88 | 81.53 | 5.75 | 109.45 | 0.95 | 30.44 | 34.21 | NA | 21.382220 | 7.061147 | NA | NA | NA | NA | NA | NA | NA |
| 0 | 21 | 86.99 | 91.33 | 93.65 | 115.00 | 123.01 | 66.84 | 10.23 | 3.48 | 11.97 | 132.15 | 90.66 | 2.58 | 68.60 | 1.04 | 30.10 | 34.00 | 48.83394 | 7.785552 | 20.442332 | NA | NA | NA | NA | NA | NA | NA |
| 0 | 16 | 114.25 | NA | NA | NA | NA | NA | 8.14 | 5.68 | 12.25 | 149.89 | NA | NA | 91.22 | 1.28 | 35.86 | 40.11 | NA | 18.489371 | 10.365972 | 1.5581016 | NA | NA | 1.159770 | NA | 0.6242926 | 0.0579575 |
| 0 | 13 | 74.80 | NA | NA | NA | NA | NA | 93.32 | 4.77 | 8.88 | 164.33 | NA | 4.33 | 84.49 | 1.07 | 23.16 | 27.76 | NA | 18.749210 | 9.068513 | NA | NA | NA | NA | NA | NA | NA |
| 1 | 17 | 84.13 | 100.23 | 69.62 | 137.84 | 91.24 | 78.76 | 28.80 | 4.22 | 14.38 | 151.67 | NA | 5.23 | 91.19 | 1.75 | 18.24 | 37.76 | 113.35402 | 20.327740 | 21.101902 | NA | NA | 1.180101 | NA | NA | 1.0568568 | NA |
| 0 | 31 | 128.68 | NA | 99.23 | 173.48 | 144.55 | 101.62 | 21.34 | 5.61 | 14.16 | 147.96 | 104.95 | 4.06 | 91.39 | 1.23 | 31.89 | 32.58 | 118.10758 | 17.200513 | 14.212995 | NA | NA | NA | 0.972422 | 1.642294 | NA | 0.5527762 |
| 0 | 20 | 110.69 | 88.69 | 72.28 | 143.15 | 193.88 | 86.83 | 25.79 | 4.50 | 13.12 | 134.66 | 68.58 | 3.56 | 81.48 | 1.11 | 28.27 | 31.64 | 225.00452 | 17.789336 | 20.238253 | NA | NA | NA | NA | NA | NA | NA |
| 0 | 16 | 83.73 | NA | NA | NA | NA | NA | 23.32 | 4.16 | 12.12 | 137.10 | 83.00 | 5.53 | 71.56 | 0.99 | 32.46 | 35.09 | 133.40626 | 20.677244 | 15.434748 | NA | NA | NA | NA | NA | NA | NA |
El pre-procesamiento es clave en el proceso de creación de un modelos de machine. Cuando queremos hacer un data frame para modelizar en aprendizaje automático, no viene todo arreglado. Incluso con una recopilación de datos buena puede que necsitemos algunos arreglos. Algunas de las causas de esto pueden ser:
Hacer todo estos procesos puede ser engorroso, pero podemos salir del paso fácilmente gracias por ejemplo de doplyr. Sin embargo ¿Cómo hacemos cuando queramos poner un modelo a grane scala y entren datos a cada momento? ¿Debemos ejecutar el pre-procesamiento una y otra vez?
La clave para resolver esto está en el paquete recipes. Esta librería trae consigo todo un conjunto de funciones que permiten un pre-procesamiento fino fino.
Para poder elaborar un flujo completo de pre-procesamiento debemos seguir los siguuientes pasos:
Hagamos un ejemplo básico: queremos una receta que sencillamnete impute los datos faltante por la media correspondiente a cada variable. Para eso vamos con ejemplo sencillo. La siguiente receta dejará los datos preparados para un modelo de aprendizaje automático, pero sin aplicarle ninguna trasnformación.
Training_datos %>%
# 🔴 1. Recipe (crea la formula generalizada del modelo)
recipe(
as.formula("Var_respuesta ~ Cognicion_01 +Cognicion_02+Escala_01+Escala_02"),
data = . ,
strata = Var_respuesta ) %>%
# 🔴 2. step
# 🟡 En este caso no aplicamos ninguina trasnformación
# STtep_log
# 🔴 3. Prep (asienta la receta y generalizala)
prep() %>%
# 🔴 4. Bake (pon la receta en producción/comprueba que tal ha funcionado)
bake(., new_data=NULL) %>%
# displayment
head(.,10) %>%
knitr::kable()| Cognicion_01 | Cognicion_02 | Escala_01 | Escala_02 | Var_respuesta |
|---|---|---|---|---|
| 66.78 | NA | 133.149284 | 18.116769 | 0 |
| 108.18 | 84.49 | 68.943373 | 12.183308 | 0 |
| 0.00 | NA | 2.898493 | 11.236520 | 0 |
| 52.88 | NA | 99.600374 | 25.026820 | 0 |
| 109.25 | 88.04 | 0.000000 | 16.798451 | 0 |
| 105.60 | 104.09 | 51.087736 | 6.731747 | 0 |
| 91.73 | NA | 160.752695 | 17.851526 | 0 |
| 91.97 | 95.35 | 15.319511 | 14.058768 | 0 |
| 69.65 | NA | 76.930164 | 22.999587 | 0 |
| 111.50 | 58.99 | 7.373436 | 8.195322 | 0 |
Training_datos %>%
# 🔴 1. Recipe (crea la formula generalizada del modelo)
recipe(
as.formula("Var_respuesta ~ Cognicion_01 +Cognicion_02+Escala_01+Escala_02"),
data = . ,
strata = Var_respuesta ) %>%
# 🔴 2. Step (haz la transofrmación que requieras)
# 🟡 Este step será el logaritmo de todas las variables predictoras númericas.
step_log(all_numeric(), -all_outcomes()) %>%
# 🔴 3. Prep (asienta la receta y generalizala)
prep() %>%
# 🔴 4. Bake (pon la receta en producción/comprueba que tal ha funcionado)
bake(., new_data=NULL) %>%
# displayment
head(.,10) %>%
knitr::kable()| Cognicion_01 | Cognicion_02 | Escala_01 | Escala_02 | Var_respuesta |
|---|---|---|---|---|
| 4.201404 | NA | 4.891471 | 2.896838 | 0 |
| 4.683796 | 4.436633 | 4.233286 | 2.500067 | 0 |
| -Inf | NA | 1.064191 | 2.419169 | 0 |
| 3.968025 | NA | 4.601166 | 3.219948 | 0 |
| 4.693639 | 4.477791 | -Inf | 2.821287 | 0 |
| 4.659658 | 4.645256 | 3.933544 | 1.906835 | 0 |
| 4.518849 | NA | 5.079867 | 2.882089 | 0 |
| 4.521462 | 4.557554 | 2.729127 | 2.643246 | 0 |
| 4.243483 | NA | 4.342898 | 3.135476 | 0 |
| 4.714025 | 4.077368 | 1.997884 | 2.103563 | 0 |
¡WALA! Datos pre-procesados y listos.
A partir de aquí en cielo es el límite ya que hay steps para todo, recomiendo leer la librería recipes y ver que tienen disponible. al final de esta sección dejaremos lista una receta como ejemplo
Hacer una imputación simple suele no ser la mejor via de hacerlo. lo ideal siempre es imputación múltiple. Es decir, una imputación que rellena los valores que faltan con números plausibles derivados de las distribuciones y relaciones entre las variables observadas en el conjunto de datos.
En tidymodels se puede hacer imputación múltiple, y valga la redundancia, de múltiples formas. Por ejemplo, podemos imputar una variable usando una regresión lineal a partir de otras. Hagamos una prueba con la variable Cognicion_02 a partir de Cognicion_01 y Escala_02
Training_datos %>%
# 🔴 1. Recipe (crea la formula generalizada del modelo)
recipe(
as.formula("Var_respuesta ~ Cognicion_01 +Cognicion_02+Escala_01+Escala_02"),
data = . ,
strata = Var_respuesta ) %>%
# 🔴 2. Step (haz la transofrmación que requieras)
# 🟡 Este step será la imputación por regresiópn lineal con otras covariables
step_impute_linear( Cognicion_02, impute_with = imp_vars(Cognicion_01, Escala_02) ) %>%
# 🔴 3. Prep (asienta la receta y generalizala)
prep() %>%
# 🔴 4. Bake (pon la receta en producción/comprueba que tal ha funcionado)
bake(., new_data=NULL) %>%
# displayment
head(.,10) %>%
knitr::kable()| Cognicion_01 | Cognicion_02 | Escala_01 | Escala_02 | Var_respuesta |
|---|---|---|---|---|
| 66.78 | 95.30378 | 133.149284 | 18.116769 | 0 |
| 108.18 | 84.49000 | 68.943373 | 12.183308 | 0 |
| 0.00 | 105.71896 | 2.898493 | 11.236520 | 0 |
| 52.88 | 97.79188 | 99.600374 | 25.026820 | 0 |
| 109.25 | 88.04000 | 0.000000 | 16.798451 | 0 |
| 105.60 | 104.09000 | 51.087736 | 6.731747 | 0 |
| 91.73 | 91.30365 | 160.752695 | 17.851526 | 0 |
| 91.97 | 95.35000 | 15.319511 | 14.058768 | 0 |
| 69.65 | 95.03224 | 76.930164 | 22.999587 | 0 |
| 111.50 | 58.99000 | 7.373436 | 8.195322 | 0 |
Esto va genial cuando conocemos la estructura de los datos y sabemos con certeza la dependencia de unas variables con otras para imputarlas fácilmente, en lugar de confiar ciegamente en algoritmos de imputación.
Como parte del pre-procesamiento, a veces puede ser necesario convertir ciertas variables a PCA para usarlas como componentes. No entraré en su uso, solo especificar la forma de hacerlo.
Training_datos %>%
# 🔴 1. Recipe (crea la formula generalizada del modelo)
recipe(
as.formula("Var_respuesta ~ Cognicion_01 + Cognicion_02 + Cognicion_03 + Escala_01 + Escala_02 + Escala_03"),
data = . ,
strata = Var_respuesta ) %>%
# 🔴 2. Step (haz la transofrmación que requieras)
# 🟡 Este step será la imputación por regresiópn lineal mcon otras covariables
step_impute_bag(all_numeric()) %>%
step_pca( Cognicion_01 , Cognicion_02 , Cognicion_03, num_comp = 2, id = "pca") %>%
# 🔴 3. Prep (asienta la receta y generalizala)
prep() %>%
# 🔴 4. Bake (pon la receta en producción/comprueba que tal ha funcionado)
tidy(id = "pca") %>%
# displayment
head(.,10) %>%
knitr::kable()| terms | value | component | id |
|---|---|---|---|
| Cognicion_01 | -0.5964794 | PC1 | pca |
| Cognicion_02 | -0.6133048 | PC1 | pca |
| Cognicion_03 | -0.5177543 | PC1 | pca |
| Cognicion_01 | 0.7421220 | PC2 | pca |
| Cognicion_02 | -0.6671338 | PC2 | pca |
| Cognicion_03 | -0.0647105 | PC2 | pca |
| Cognicion_01 | 0.3057241 | PC3 | pca |
| Cognicion_02 | 0.4228354 | PC3 | pca |
| Cognicion_03 | -0.8530785 | PC3 | pca |
ADASYN, o Adaptive Synthetic Sampling, es un algoritmo bastante utilizado cuando tenemos conjuntos de datos desequilibrados.
Los conjuntos de datos desequilibrados se producen cuando la distribución de clases en los datos es desigual, lo que dificulta el entrenamiento de modelos que pueden estar sesgados hacia la clase mayoritaria.
# 🔴 Libreria que permite instalar algoritmos e sobremuestreo en tidymodels.
pacman::p_load(themis)
Training_datos_adasyn <- Training_datos %>%
# 1. Recipe (crea la formula generalizada del modelo)
recipe(
as.formula("Var_respuesta ~ Cognicion_01 +Cognicion_02+Escala_01+Escala_02"),
data = . ,
strata = Var_respuesta ) %>%
# 2. Step (haz la transofrmación que requieras)
# step de imputación por arboles manteniendo la estrucura de todas las covariables
step_impute_bag(all_numeric()) %>%
# 🔴 step para sobremuestrear la el nivel menor en la variable respuesta.
step_adasyn(Var_respuesta, over_ratio = 1, neighbors = 10) %>%
# 3. Prep (asienta la receta y generalizala)
prep() %>%
# 4. Bake (pon la receta en producción/comprueba que tal ha funcionado)
bake(., new_data=NULL)| Cognicion_01 | Cognicion_02 | Escala_01 | Escala_02 | Var_respuesta |
|---|---|---|---|---|
| 66.78 | 102.05505 | 133.149284 | 18.116769 | 0 |
| 108.18 | 84.49000 | 68.943373 | 12.183308 | 0 |
| 0.00 | 99.58442 | 2.898493 | 11.236520 | 0 |
| 52.88 | 97.03067 | 99.600374 | 25.026820 | 0 |
| 109.25 | 88.04000 | 0.000000 | 16.798451 | 0 |
| 105.60 | 104.09000 | 51.087736 | 6.731747 | 0 |
| 91.73 | 89.86937 | 160.752695 | 17.851526 | 0 |
| 91.97 | 95.35000 | 15.319511 | 14.058768 | 0 |
| 69.65 | 100.57019 | 76.930164 | 22.999587 | 0 |
| 111.50 | 58.99000 | 7.373436 | 8.195322 | 0 |
ggpubr::ggarrange(
Training_datos %>%
group_by(Var_respuesta) %>%
summarise(recuento = n()) %>%
ungroup() %>%
mutate(prop_var_respuesta = recuento / sum(recuento) ) %>%
ggplot(aes(x = Var_respuesta, y = recuento, fill = Var_respuesta)) +
geom_bar(stat = "identity") +
geom_text(
aes(label = scales::percent(prop_var_respuesta)),
position = position_stack(vjust = 0.5), size = 4) +
scale_fill_manual( values=c("azure4", "azure3"))
,
Training_datos_adasyn %>%
group_by(Var_respuesta) %>%
summarise(recuento = n()) %>%
mutate(prop_var_respuesta = recuento / sum(recuento) ) %>%
ungroup() %>%
ggplot(aes(x = Var_respuesta, y = recuento, fill = Var_respuesta)) +
geom_bar(stat = "identity") +
geom_text(
aes(label = scales::percent(prop_var_respuesta)),
position = position_stack(vjust = 0.5), size = 4) +
scale_fill_manual( values=c("azure4", "azure3")),
common.legend = T, legend = 'bottom'
,
labels = c('Normal', 'Sobre Muestreo\n Adasyn'), label.x = 0.1
) %>%
plot()Por el momento dejaré una receta lista para continuar con el taller.
# Creamos la Fórmula para la receta a partir
# de nombres de las variables y un poco de magia con paste
Receta_modelo_1_formula <- Training_datos %>%
dplyr::select(
dplyr::matches('Cog'), dplyr::matches('escala')) %>%
names() %>%
paste(., collapse = ' + ') %>%
paste('Var_respuesta ~' ,.) %>%
as.formula()
# 🔴 1.) Iniciamos la receta
Receta_modelo_1 <- recipe(
formula = Receta_modelo_1_formula,
data = Training_datos ,
strata = Var_respuesta) %>%
# 🔴 2.) Steps
# 🟡 2.1) Eliminamos las variables que contengan más de un 20% de datos perdidos
step_filter_missing(all_predictors(), threshold = 0.1) %>%
# 🟡 2.2) Imputamos las variables númericas con un algoritmo de bagged trees
step_impute_bag(all_numeric(),-all_outcomes()) %>%
# 2.3) Normalizamos los datos (restamos la media)
step_normalize(all_numeric(),-all_outcomes() ) %>%
# 2.4) Escalamos los datos (reducimos a escala entre 0 y 1)
step_scale(all_numeric(),-all_outcomes()) %>%
# 2.5) Convertimos a Dummy las variables factor (no hará nada por que nop tenemos variables factor)
step_dummy(all_nominal(), -all_outcomes(), one_hot = T) %>%
# 2.6) Sobre muestreamos la variable para equilibrar grupos
step_adasyn(Var_respuesta, over_ratio = 1, neighbors = 10) %>%
# 3.) Preparamos la receta
prep()| Cognicion_01 | Escala_02 | Escala_03 | Var_respuesta |
|---|---|---|---|
| -0.6426710 | -0.0464912 | 1.5341196 | 0 |
| 0.7604689 | -0.7638606 | -0.3800329 | 0 |
| -2.9059967 | -0.8783295 | -1.4610980 | 0 |
| -1.1137735 | 0.7889501 | -0.7053657 | 0 |
| 0.7967336 | -0.2058789 | -1.3439820 | 0 |
| 0.6730268 | -1.4229671 | -0.9510280 | 0 |
| 0.2029411 | -0.0785598 | -0.0393877 | 0 |
| 0.2110752 | -0.5371131 | -0.1631135 | 0 |
| -0.5454002 | 0.5438529 | -0.4982340 | 0 |
| 0.8729912 | -1.2460174 | -1.4502475 | 0 |
Vale, dejamos atrás la parte del pre-procesamiento y empezamos la parte más machine leaninrg. Modelos. Modelos a la carrera.
La función del paquete parsnip dentro de tidymodels es propoconar una interfaz universal para diferentes métodos de modelado en R.
El ejemplo del principio del taller está basado en esta sección, parnsip es el encargado de que solo con cambiar unas pocas líneas puedas cambiar todo un modelo dentro de un flujo de aprendizaje automático.
A fecha de este taller, Tidymodels cuenta con más de 155 modelos, incluyendo supervivencia y series de tiempo. Se pueden encontrar en: https://www.tidymodels.org/find/parsnip/
Supongamos que queremos empezar modelar, queremos usar una regresión logística, un clásico en este campo.Pero dentro de R hay varias librarías que usan la regresión logística, cada una con sus varianciones y reglas. Veamos que tenemos para esto:
## # A tibble: 8 × 2
## engine mode
## <chr> <chr>
## 1 glm classification
## 2 glmnet classification
## 3 LiblineaR classification
## 4 spark classification
## 5 keras classification
## 6 stan classification
## 7 brulee classification
## 8 h2o classification
logistic_reg(
engine = "glmnet",
#Hyper-parametros
penalty = NULL, mixture = NULL, mode = 'classification')## Logistic Regression Model Specification (classification)
##
## Computational engine: glmnet
logistic_reg(
engine = "glmnet",
#Hyper-parametros
penalty = 0, mixture = 0, mode = 'classification')## Logistic Regression Model Specification (classification)
##
## Main Arguments:
## penalty = 0
## mixture = 0
##
## Computational engine: glmnet
logistic_reg(
engine = "glmnet",
#Hyper-parametros
penalty = tune(), mixture = tune(), mode = 'classification')## Logistic Regression Model Specification (classification)
##
## Main Arguments:
## penalty = tune()
## mixture = tune()
##
## Computational engine: glmnet
Model_RandomForest <-
# 🔴 Especificamos el modelo que queremos en este caso un random forest
rand_forest(
# 🟡 Ajustamos los hiper-parametros
mtry = tune(), trees = tune(), min_n = tune()) %>%
# 🔴 Ponemos el motor, es decir la librería por la queremos que se ejecute el modelo
set_engine("ranger", importance = "impurity") %>%
# 🔴 Ajustamos el modo del modelo, es decri si queremos regresion o clasifiacion
set_mode("classification")El paquete de Workflows permite combinar una receta y una especificación de modelo en un único objeto. Esto hace que todo el proceso de modelado, desde el preprocesamiento de datos hasta el ajuste del modelo, esté en un mismo punto
El paquete yardstick, en concreto, se utiliza para las métricas de evaluación de modelos. Además
# 🔴 Seleccionamos las metricas que queremos para nuestro modelo
Modelo_Metricas <- metric_set(accuracy, j_index, precision, sensitivity, specificity, roc_auc, f_meas,recall, mcc)También decir que es es posible crear tu propia métrica. Pero aún no he indagado en esta metodología. Se puede encontrar más en: https://yardstick.tidymodels.org/articles/custom-metrics.html
Los hiper-parámetros son una pieza crucial en los modelos de aprendizaje automático. Estos son parámetros de configuración externos que no pueden aprenderse a partir de los datos, pero que influyen significativamente en el rendimiento del modelo. Un ejemplo sería el número de árboles de un random forest por ejemplo. A diferencia de los parámetros del modelo, que se aprenden durante el entrenamiento, los híper-parámetros se establecen antes de que comience el entrenamiento de este.
La selección de los hiperparámetros adecuados es esencial para conseguir un modelo de aprendizaje automático generalizable y con buen rendimiento y es un proceso delicado. Los híper-parámetros desempeñan un papel crucial a la hora de controlar el equilibrio entre la sobreadaptación y la inadaptación (oiverfitting y underfitting). La sobreadaptación se produce cuando un modelo funciona bien con los datos de entrenamiento pero no consigue generalizarse a datos nuevos que no se han visto. Por otro lado, el infraajuste se produce cuando el modelo es demasiado simple y no puede captar los patrones subyacentes en los datos. Ajustar los hiperparámetros puede ayudar a encontrar el nivel adecuado de complejidad del modelo, haciendo un buen equilibrio entre la varianza del conjunto de entrenamiento y el de evaluación.
Los paquetes Tune y Dials son los que nos permiten hacer esta selección y control. Aquí se crea todo un subflujo de trabajo para el ajuste de híper-parámetros.
Grid regular hará las combinaciones acorde al rango que le hayamos dado. Por defecto hace la cmbinación entre el mínimo, el medio y el máximo de cada valor.
## # A tibble: 27 × 3
## mtry min_n trees
## <int> <int> <int>
## 1 2 2 2
## 2 6 2 2
## 3 10 2 2
## 4 2 6 2
## 5 6 6 2
## 6 10 6 2
## 7 2 10 2
## 8 6 10 2
## 9 10 10 2
## 10 2 2 6
## # ℹ 17 more rows
Podemos ampliar esto con el parámetro levels. Ahora serrá con el mínimo, el quantile 33, el cuantil 66 y el máximo.
grid_regular(
mtry(range = c(2, 10)),
min_n(range = c(2, 10)),
trees(range = c(2, 10)),
levels = 4 )## # A tibble: 64 × 3
## mtry min_n trees
## <int> <int> <int>
## 1 2 2 2
## 2 4 2 2
## 3 7 2 2
## 4 10 2 2
## 5 2 4 2
## 6 4 4 2
## 7 7 4 2
## 8 10 4 2
## 9 2 7 2
## 10 4 7 2
## # ℹ 54 more rows
grid_max_entropy(
mtry(range = c(1, 4)),
min_n(range = c(10, 30)),
trees(range = c(1, 1000)),
size = 10
)## # A tibble: 10 × 3
## mtry min_n trees
## <int> <int> <int>
## 1 1 25 20
## 2 2 11 610
## 3 4 12 2
## 4 2 20 936
## 5 4 12 733
## 6 4 30 868
## 7 1 20 676
## 8 2 19 307
## 9 2 29 661
## 10 4 26 190
WF_hiper_parametros <-
# 🔴 Este proceso se hace con la función tune grid 🔴
tune_grid(
# Ponemos la receta nuestro modelo, que incuye el preprocesamiento y la fórmula
object = WF_Receta_modelo_1_Random_forest,
# 🟡 Ponemos los Folds del conjunto de entrenamiento,
# 🟡 Dentro del entrenamiento hará un sub entrenamiento y una sub validación
resamples = Folds_training,
#Ponemos un grid para que pruebe con diferentes híper-parámetros, ya con un rango pre-definido
grid = grid_max_entropy(
mtry(range = c(1, 4)),
min_n(range = c(10, 30)),
trees(range = c(1, 1000)),
size = 10),
# Ponemos las métricas que hemos definido anteriormente
metrics = Modelo_Metricas,
# Con esta opción podemos guardar las predicciones
control = control_grid( save_pred = T)
)WF_hiperparametros_coleccion <-
WF_hiper_parametros %>%
collect_metrics() # 🔴 Esta es la funcion que retorna las métricas, importante tenerla en cuenta
WF_hiperparametros_coleccion %>%
head(.,20) %>%
knitr::kable()| mtry | trees | min_n | .metric | .estimator | mean | n | std_err | .config |
|---|---|---|---|---|---|---|---|---|
| 2 | 886 | 10 | accuracy | binary | 0.6711538 | 5 | 0.0523962 | Preprocessor1_Model01 |
| 2 | 886 | 10 | f_meas | binary | 0.7730576 | 5 | 0.0379652 | Preprocessor1_Model01 |
| 2 | 886 | 10 | j_index | binary | 0.2332661 | 5 | 0.1447212 | Preprocessor1_Model01 |
| 2 | 886 | 10 | mcc | binary | 0.2023982 | 5 | 0.1231124 | Preprocessor1_Model01 |
| 2 | 886 | 10 | precision | binary | 0.8541935 | 5 | 0.0368454 | Preprocessor1_Model01 |
| 2 | 886 | 10 | recall | binary | 0.7082661 | 5 | 0.0439273 | Preprocessor1_Model01 |
| 2 | 886 | 10 | roc_auc | binary | 0.7237021 | 5 | 0.0503109 | Preprocessor1_Model01 |
| 2 | 886 | 10 | sensitivity | binary | 0.7082661 | 5 | 0.0439273 | Preprocessor1_Model01 |
| 2 | 886 | 10 | specificity | binary | 0.5250000 | 5 | 0.1145644 | Preprocessor1_Model01 |
| 4 | 123 | 28 | accuracy | binary | 0.6610256 | 5 | 0.0440172 | Preprocessor1_Model02 |
| 4 | 123 | 28 | f_meas | binary | 0.7638920 | 5 | 0.0305758 | Preprocessor1_Model02 |
| 4 | 123 | 28 | j_index | binary | 0.2393145 | 5 | 0.1364279 | Preprocessor1_Model02 |
| 4 | 123 | 28 | mcc | binary | 0.2009226 | 5 | 0.1160916 | Preprocessor1_Model02 |
| 4 | 123 | 28 | precision | binary | 0.8594062 | 5 | 0.0349711 | Preprocessor1_Model02 |
| 4 | 123 | 28 | recall | binary | 0.6893145 | 5 | 0.0333552 | Preprocessor1_Model02 |
| 4 | 123 | 28 | roc_auc | binary | 0.7299773 | 5 | 0.0469362 | Preprocessor1_Model02 |
| 4 | 123 | 28 | sensitivity | binary | 0.6893145 | 5 | 0.0333552 | Preprocessor1_Model02 |
| 4 | 123 | 28 | specificity | binary | 0.5500000 | 5 | 0.1159202 | Preprocessor1_Model02 |
| 2 | 81 | 16 | accuracy | binary | 0.6761538 | 5 | 0.0550452 | Preprocessor1_Model03 |
| 2 | 81 | 16 | f_meas | binary | 0.7766291 | 5 | 0.0400532 | Preprocessor1_Model03 |
WF_hiper_parametros %>%
collect_metrics() %>%
pivot_longer(cols = c('mtry','trees', 'min_n'), names_to = 'tipo_parametro', values_to = 'valor_parametro') %>%
ggplot(aes(valor_parametro,mean, color=mean)) +
geom_point() +
facet_grid(.metric~tipo_parametro, scales = 'free' ) +
scale_color_viridis_b()WF_hiperparametros_mejor <- WF_hiper_parametros %>% show_best("roc_auc",1)
WF_hiperparametros_mejor %>%
knitr::kable()| mtry | trees | min_n | .metric | .estimator | mean | n | std_err | .config |
|---|---|---|---|---|---|---|---|---|
| 2 | 240 | 29 | roc_auc | binary | 0.734816 | 5 | 0.0384686 | Preprocessor1_Model06 |
Ya hemos hecho todos los procesos previos. Todo lo que a un modelo de machine learning le puede hacer falta. Finalmente, una vez que lo tenemos todo hilado, hay que proceder a entrenar el modelo y dejarlo listo.
Modelo_fitted <- WF_Receta_modelo_1_Random_forest %>%
# 🔴 Finalizamos el flujo de trabajo con la mejor combinación de Hiper paraemtros encontrada
finalize_workflow(WF_hiperparametros_mejor) %>%
# 🔴 finalmente ajustamos con fit
fit(data= Training_datos)Con fit ya por fin tenemos nuestro modelo, ya podemos ponernos a hacer lo que hace el aprendizaje automático: predecir. vamos a crear una observación nueva y a ver como la trabaja.
Nuevo_caso <- tibble(
'Cognicion_01' = 123,
'Cognicion_02' = 102,
'Cognicion_03' = 89,
'Cognicion_04' = 100,
'Cognicion_05' = 179,
'Escala_01' = 12,
'Escala_02' = 22,
'Escala_03' = 20
)tibble(
# 🔴 La función predict de la librería stats tiene sinergias con Tidymodels, neecsitra de estos argumentos:
# - 🟡 El modelo ya ajustado (model_fitted)
# - 🟡 Un tibble/dataframe con la nueva prediccion (debe contener las mismas variables)
# - 🟡 Un tipo de prediccion en los casos de clasificación, si queremos las probabiliadades de cada categoria o la categoría predicha directamente
predict(Modelo_fitted, Nuevo_caso, type = "prob" ),
predict(Modelo_fitted, Nuevo_caso, type = "class")) %>%
knitr::kable()| .pred_0 | .pred_1 | .pred_class |
|---|---|---|
| 0.9607399 | 0.0392601 | 0 |
tibble(predict(Modelo_fitted, Testing_datos, type = "prob" ),predict(Modelo_fitted, Testing_datos, type = "class")) %>%
head(.,10) %>%
knitr::kable()| .pred_0 | .pred_1 | .pred_class |
|---|---|---|
| 0.7925329 | 0.2074671 | 0 |
| 0.6064286 | 0.3935714 | 0 |
| 0.8727562 | 0.1272438 | 0 |
| 0.4408074 | 0.5591926 | 1 |
| 0.8695524 | 0.1304476 | 0 |
| 0.9175403 | 0.0824597 | 0 |
| 0.0714041 | 0.9285959 | 1 |
| 0.9023667 | 0.0976333 | 0 |
| 0.9603678 | 0.0396322 | 0 |
| 0.7704542 | 0.2295458 | 0 |
Crear las predcciones y todo es un sistema engorroso. Además de que se pierde mucha informaciíon por el camino, ya que todo sigue estando en objetos dispersos. Para ello los desarrolladores crearon last_fit.
La función last_fit permite tener todo unificado en un solo tibble, al tener que ajustarse con un objeto split, esta recibe directamente los datos de entrenamiento y los datos de evaulación (revibe datos de training y testing) y entonces crea sus métricas y sus predicciones.
Con lsat_fit contamos con una serie de funciones para extraer todo lo que necesitemos de ellas: - extract_fit_engine() : permite extraer el modelo crudo como si no hubiera sido ajustado con tidymodels, de modo que si es un modelo explicativo puede ser tratado con normalidad - collect_metrics() : extraer un tibble con las métricas del modelo (puestas con metric set) - collect_predictions() : extrae un tibble con las predicciones de los datos de evaluacion. (Por defecto el umbral de rendimiento en esta función es 0.5, y estoy buscando alguna forma de cambiar esto, pero se puede manipular el mismo tibble como cualquier otro, de modo que no es un drama.) - extract_workflow() : Permite extraer todo el flujo del modelo, de incluyendo la receta de modo que podemos hacer predicciones:
Modelo_1_final <- WF_Receta_modelo_1_Random_forest %>%
# 🔴 Finalizamos el flujo de trabajo con la mejor combinación de Hiper paraemtros encontrada
finalize_workflow(WF_hiperparametros_mejor) %>%
# 🔴 Ajustamos con last_fit, poniendo el split de datos inicial y las metricas para las predicciones
last_fit(Split_datos_strat, metrics = Modelo_Metricas)
Modelo_1_final## # Resampling results
## # Manual resampling
## # A tibble: 1 × 6
## splits id .metrics .notes .predictions .workflow
## <list> <chr> <list> <list> <list> <list>
## 1 <split [198/87]> train/test split <tibble> <tibble> <tibble> <workflow>
Modelo_1_final %>%
# 🔴 dentro del flujo de trabajo están la receta y el modelo ya ajustado
extract_workflow() %>%
# de modo que lo podemos usar para hacer predicciones
predict(Testing_datos, type = "prob" )## # A tibble: 87 × 2
## .pred_0 .pred_1
## <dbl> <dbl>
## 1 0.793 0.207
## 2 0.606 0.394
## 3 0.873 0.127
## 4 0.441 0.559
## 5 0.870 0.130
## 6 0.918 0.0825
## 7 0.0714 0.929
## 8 0.902 0.0976
## 9 0.960 0.0396
## 10 0.770 0.230
## # ℹ 77 more rows
A partir de aquí podemos manipular el objeto fit de múltiples maneras múltiples
Modelo_1_final %>%
# 🔴 Con esta funcion podemos recopilar todas las metricas de un modelo
collect_metrics() %>%
select(-.estimator,-.config) %>%
mutate(.estimate= round(.estimate,2)) %>%
knitr::kable()| .metric | .estimate |
|---|---|
| accuracy | 0.71 |
| j_index | 0.39 |
| precision | 0.89 |
| sensitivity | 0.72 |
| specificity | 0.67 |
| f_meas | 0.80 |
| recall | 0.72 |
| mcc | 0.33 |
| roc_auc | 0.74 |
options(yardstick.event_first = FALSE)
Modelo_1_final %>% collect_predictions() %>%
# 🔴 Esta funcion permite calcular rapidamente la curva roc de un modelo ya ajustado
roc_curve(truth = Var_respuesta, .pred_0 ) %>%
ggplot(aes( x = 1- specificity ,y= sensitivity ) ) +
geom_path() +
geom_abline(intercept=0, slope=1, linetype=3) +
labs(title= 'Cruva roc del modelo_1_final')Advertencia: en la función de curva roc, ved que puesto prediccioón del 0 en vez del 1. Esto es un error interno en la clasificación dentro del del paquete yardstick. Hay que tener cuidado sobre si estas funciones utilizan el primer / segundo / x nivel de su factor como el “evento”. Por defecto, yardstick elige el primer nivel de verdad como “evento” al calcular la curva roc. Podeíos encontrar más de este suceso en el github de los desarrolladores https://github.com/tidymodels/yardstick/issues/94
El concepto de rendimiento del umbral (threshold performance) de es crucial a la hora de evaluar y ajustar modelos de clasificación binaria.
El umbral es el valor de probabilidad por encima del cual un resultado predicho se considera de clase positiva. Por defecto, muchos modelos utilizan un umbral de 0,5. Sin embargo, en casi todos los casos, el umbral óptimo no estará otro punto. Este “óptimo lo tendremo que deifnir que en función de los objetivos y requisitos específicos del problema.
Ajustar el umbral permite elegir entre precisión y la recall. Bajar el umbral aumenta la sensibilidad (recall) pero disminuye la precisión, y viceversa. Básicamente, el umbral óptimo depende de la importancia relativa de los falsos positivos y los falsos negativos en una aplicación concreta.
En tidymodels está la librería probably con la función threshold_perf, la cual nos permite, dadas las predicciones del modelo ver como se desenvolupa este umbral.
pacman::p_load(probably)
Modelo_1_final_Threshold_performance <-
# Llamamos al modelo ya finallizado con "last_fit"
Modelo_1_final %>%
# Coleccionamos sus predicciones
collect_predictions() %>%
# 🔴 Ejecutamos el redimeinto de umbral con la siguiente función
threshold_perf(
# Le decimos que variable es la verdadera respuesta
truth = Var_respuesta,
# le decimos que vairbale es la predicción (aviso: debido a un bugg en la version )
.pred_1,
# ajustamos un rango en el qeu se evaluará el umbral en cada punto
thresholds = seq(0.5,1,0.01) )
# displayment
Modelo_1_final_Threshold_performance %>%
head(.,10) %>%
mutate(.estimate= round(.estimate,2)) %>%
knitr::kable()| .threshold | .metric | .estimator | .estimate |
|---|---|---|---|
| 0.50 | sensitivity | binary | 0.28 |
| 0.51 | sensitivity | binary | 0.28 |
| 0.52 | sensitivity | binary | 0.28 |
| 0.53 | sensitivity | binary | 0.28 |
| 0.54 | sensitivity | binary | 0.26 |
| 0.55 | sensitivity | binary | 0.26 |
| 0.56 | sensitivity | binary | 0.25 |
| 0.57 | sensitivity | binary | 0.25 |
| 0.58 | sensitivity | binary | 0.25 |
| 0.59 | sensitivity | binary | 0.25 |
Modelo_1_final_Threshold_performance %>%
ggplot(aes(x=.threshold,y=.estimate ,color=.metric)) +
geom_line(size=1.1)+
scale_color_manual(values=c('purple','blue','pink'))+
geom_vline(xintercept =
Modelo_1_final_Threshold_performance %>%
filter(.metric=='j_index') %>%
arrange(desc(.estimate)) %>%
slice(1) %>%
pull(.threshold), linetype = 'dashed', size=1.5)Bien, ya hemos aprendido a crear todo un flujo entero de machine learning con Tidymodels, desde la base hasta su final. Muy útil todo. pero calibrar modelos de 1 en 1 es un auténtico petardazo. No me acaba de ser útil este taller. ¿Cómo me va a ser útil si tengo que porbar unas 20 fórmulas con datos?
Lo entiendo perfectamente, a mi también me pasó. Los médicos me preguntaron que era mejor para tener en cuenta, si las escalas, las congniciones la bioquimica o la genetica a la hora de diagnosticar al paciente.
Comprovemoslo.
Receta_modelo_Escalas_formula <- Training_datos %>%
dplyr::select(dplyr::matches('Esc')) %>%
names() %>%
paste(., collapse = ' + ') %>%
paste('Var_respuesta ~' ,.) %>%
as.formula()
# 🔴 1.) Iniciamos la receta
Receta_modelo_Escalas <- recipe(
formula = Receta_modelo_Escalas_formula,
data = Training_datos ,
strata = Var_respuesta) %>%
# 🔴 2.) Steps
# 🟡 2.1) Eliminamos las variables que contengan más de un 25% de datos perdidos
step_filter_missing(all_predictors(), threshold = 0.25) %>%
# 🟡 2.2) Imputamos las variables númericas con un algoritmo de bagged trees
step_impute_bag(all_numeric(),-all_outcomes()) %>%
# 2.3) Normalizamos los datos (restamos la media)
step_normalize(all_numeric(),-all_outcomes() ) %>%
# 2.4) Escalamos los datos (reducimos a escala entre 0 y 1)
step_scale(all_numeric(),-all_outcomes()) %>%
# 2.5) Convertimos a Dummy las variables factor (no hará nada por que nop tenemos variables factor)
step_dummy(all_nominal(), -all_outcomes(), one_hot = T) %>%
# 2.6) Sobre muestreamos la variable para equilibrar grupos
step_adasyn(Var_respuesta, over_ratio = 1, neighbors = 10) %>%
# 3.) Preparamos la receta
prep()
Receta_modelo_Cognicion_formula <- Training_datos %>%
dplyr::select(
dplyr::matches('Cog')) %>%
names() %>%
paste(., collapse = ' + ') %>%
paste('Var_respuesta ~' ,.) %>%
as.formula()
# 🔴 1.) Iniciamos la receta
Receta_modelo_Cognicion <- recipe(
formula = Receta_modelo_Cognicion_formula,
data = Training_datos ,
strata = Var_respuesta) %>%
# 🔴 2.) Steps
# 🟡 2.1) Eliminamos las variables que contengan más de un 25% de datos perdidos
step_filter_missing(all_predictors(), threshold = 0.25) %>%
# 🟡 2.2) Imputamos las variables númericas con un algoritmo de bagged trees
step_impute_mean(all_numeric(),-all_outcomes()) %>%
# 2.3) Normalizamos los datos (restamos la media)
step_normalize(all_numeric(),-all_outcomes() ) %>%
# 2.4) Escalamos los datos (reducimos a escala entre 0 y 1)
step_scale(all_numeric(),-all_outcomes()) %>%
# 2.5) Convertimos a Dummy las variables factor (no hará nada por que nop tenemos variables factor)
step_dummy(all_nominal(), -all_outcomes(), one_hot = T) %>%
# 2.6) Sobre muestreamos la variable para equilibrar grupos
step_adasyn(Var_respuesta, over_ratio = 1, neighbors = 10) %>%
# 3.) Preparamos la receta
prep()
bake(Receta_modelo_Cognicion, new_data = NULL)## # A tibble: 316 × 4
## Cognicion_01 Cognicion_03 Cognicion_04 Var_respuesta
## <dbl> <dbl> <dbl> <fct>
## 1 -0.643 0 0 0
## 2 0.760 0.313 1.13 0
## 3 -2.91 0 0 0
## 4 -1.11 0 0 0
## 5 0.797 1.64 0.259 0
## 6 0.673 0.570 0.393 0
## 7 0.203 0 0 0
## 8 0.211 0.203 0.886 0
## 9 -0.545 0 0 0
## 10 0.873 1.50 1.09 0
## # ℹ 306 more rows
Receta_modelo_Bioquimicas_formula <- Training_datos %>%
dplyr::select(
dplyr::matches('Bio')) %>%
names() %>%
paste(., collapse = ' + ') %>%
paste('Var_respuesta ~' ,.) %>%
as.formula()
# 🔴 1.) Iniciamos la receta
Receta_modelo_Bioquimica <- recipe(
formula = Receta_modelo_Bioquimicas_formula,
data = Training_datos ,
strata = Var_respuesta) %>%
# 🔴 2.) Steps
# 🟡 2.1) Eliminamos las variables que contengan más de un 25% de datos perdidos
step_filter_missing(all_predictors(), threshold = 0.25) %>%
# 🟡 2.2) Imputamos las variables númericas con un algoritmo de bagged trees
step_impute_mean(all_numeric()) %>%
# 2.3) Normalizamos los datos (restamos la media)
step_normalize(all_numeric(),-all_outcomes() ) %>%
# 2.4) Escalamos los datos (reducimos a escala entre 0 y 1)
step_scale(all_numeric(),-all_outcomes()) %>%
# 2.5) Convertimos a Dummy las variables factor (no hará nada por que nop tenemos variables factor)
step_dummy(all_nominal(), -all_outcomes(), one_hot = T) %>%
# 2.6) Sobre muestreamos la variable para equilibrar grupos
step_adasyn(Var_respuesta, over_ratio = 1, neighbors = 10) %>%
# 3.) Preparamos la receta
prep()
Receta_modelo_Popurri <- Training_datos %>%
dplyr::select(
Cognicion_01, Cognicion_02, Escala_01,Escala_02, Bioquimica_01,Bioquimica_02, Genetica_01,Genetica_02) %>%
names() %>%
paste(., collapse = ' + ') %>%
paste('Var_respuesta ~' ,.) %>%
as.formula()
# 🔴 1.) Iniciamos la receta
Receta_modelo_Popurri <- recipe(
formula = Receta_modelo_Popurri,
data = Training_datos ,
strata = Var_respuesta) %>%
# 🔴 2.) Steps
# 🟡 2.1) Eliminamos las variables que contengan más de un 5% de datos perdidos
step_filter_missing(all_predictors(), threshold = 0.25) %>%
# 🟡 2.2) Imputamos las variables númericas con un algoritmo de bagged trees
step_impute_bag(all_numeric(),-all_outcomes()) %>%
# 2.3) Normalizamos los datos (restamos la media)
step_normalize(all_numeric(),-all_outcomes() ) %>%
# 2.4) Escalamos los datos (reducimos a escala entre 0 y 1)
step_scale(all_numeric(),-all_outcomes()) %>%
# 2.5) Convertimos a Dummy las variables factor (no hará nada por que nop tenemos variables factor)
step_dummy(all_nominal(), -all_outcomes(), one_hot = T) %>%
# 2.6) Sobre muestreamos la variable para equilibrar grupos
step_adasyn(Var_respuesta, over_ratio = 1, neighbors = 10) %>%
# 3.) Preparamos la receta
prep()Ahora, esto son muchas recetas y mucho texto. Será dificil combinarlo todo, ¿no?.¡Que va! Workflow_set ya está pensado para esto. Funciona como el crossing de la libreia purrr o como el expand.grid de la librería base, crear workflows con todas las combinaciones.
Receta_list <- list(
'Modelo_Escalas' = Receta_modelo_Escalas,
'Modelo_Cognicion' = Receta_modelo_Cognicion,
'Modelo_Bioquimica' = Receta_modelo_Bioquimica,
'Modelo_Popurri' = Receta_modelo_Popurri
)
Model_list <- list(
RandomForest = Model_RandomForest
# ,XGBoost = Model_XGBoost
# ,Bagged_Trees = Model_Bagged_Tree
)
Wokflows_set <- workflow_set(
preproc = Receta_list,
models = Model_list,
# 🔴 La opcion Cross es para que haga todos los cruces de modelos con recetas
cross = T)
Wokflows_set## # A workflow set/tibble: 4 × 4
## wflow_id info option result
## <chr> <list> <list> <list>
## 1 Modelo_Escalas_RandomForest <tibble [1 × 4]> <opts[0]> <list [0]>
## 2 Modelo_Cognicion_RandomForest <tibble [1 × 4]> <opts[0]> <list [0]>
## 3 Modelo_Bioquimica_RandomForest <tibble [1 × 4]> <opts[0]> <list [0]>
## 4 Modelo_Popurri_RandomForest <tibble [1 × 4]> <opts[0]> <list [0]>
Una vez hecho el workfñow_set, se pone a trabajr con workflow_map.
Wokflows_set_map <-
Wokflows_set %>%
# 🔴 workflow maps es una funcion que permite ajustar multiples flujos de trabajo
workflow_map(
resamples = Folds_training,
fn = "tune_grid",
grid = grid_max_entropy(
mtry(range = c(1, 4)),
min_n(range = c(10, 30)),
trees(range = c(1, 1000)),
size = 10),
# verbose = TRUE,
metrics = Modelo_Metricas,
control = control_grid( save_pred = T),
# 🔴 para garantizare replicabildiad, podemos fijar la semilla en el proceso
seed = 2465)
Wokflows_set_map## # A workflow set/tibble: 4 × 4
## wflow_id info option result
## <chr> <list> <list> <list>
## 1 Modelo_Escalas_RandomForest <tibble [1 × 4]> <opts[4]> <tune[+]>
## 2 Modelo_Cognicion_RandomForest <tibble [1 × 4]> <opts[4]> <tune[+]>
## 3 Modelo_Bioquimica_RandomForest <tibble [1 × 4]> <opts[4]> <tune[+]>
## 4 Modelo_Popurri_RandomForest <tibble [1 × 4]> <opts[4]> <tune[+]>
AVISO, workflow_map puede paralelizarse, pero aún estoy descubriendo como integrarlo de forma efectiva.
Finalmnete, podemos calibrar cada modelo individualmente y tener todo finiquitado. Ese proceso sería más eficiente con una funcion. desde luego. Pero se me acaba el tiempo para dedicar al taller :(
# Modelo Escalas ----
Modelo_Escalas_RandomForest_hyperparametros <- Wokflows_set_map %>%
# 🔴 Extraemos el flujo correspondiente al modelo que queremos
extract_workflow_set_result('Modelo_Escalas_RandomForest') %>%
# 🔴 podemos ver una tabla con los mejores resultados, acorde a una metrica
show_best(metric = 'sensitivity', n=50) %>%
filter(row_number()==1) %>%
select(mtry,trees,min_n,.config )
Modelo_Escalas_RandomForest_final <- Wokflows_set_map %>%
# 🔴 Extraemos el flujo correspondiente al modelo que queremos
extract_workflow('Modelo_Escalas_RandomForest') %>%
finalize_workflow(Modelo_Escalas_RandomForest_hyperparametros) %>%
last_fit(Split_datos_strat, metrics = Modelo_Metricas )
# Modelo_Cognicion ----
Modelo_Cognicion_RandomForest_hyperparametros <- Wokflows_set_map %>%
extract_workflow_set_result('Modelo_Cognicion_RandomForest') %>%
show_best(metric = 'sensitivity', n=50) %>%
filter(row_number()==1) %>%
select(mtry,trees,min_n,.config )
Modelo_Cognicion_RandomForest_final <- Wokflows_set_map %>%
extract_workflow('Modelo_Cognicion_RandomForest') %>%
finalize_workflow(Modelo_Cognicion_RandomForest_hyperparametros) %>%
last_fit(Split_datos_strat, metrics = Modelo_Metricas )
# Modelo_Bioquimico ----
Modelo_Bioquimica_RandomForest_hyperparametros <- Wokflows_set_map %>%
extract_workflow_set_result('Modelo_Bioquimica_RandomForest') %>%
show_best(metric = 'sensitivity', n=50) %>%
filter(row_number()==1) %>%
select(mtry,trees,min_n,.config )
Modelo_Bioquimica_RandomForest_final <- Wokflows_set_map %>%
extract_workflow('Modelo_Bioquimica_RandomForest') %>%
finalize_workflow(Modelo_Bioquimica_RandomForest_hyperparametros) %>%
last_fit(Split_datos_strat, metrics = Modelo_Metricas )
# Modelo Popurri ----
Modelo_Popurri_RandomForest_hyperparametros <- Wokflows_set_map %>%
extract_workflow_set_result('Modelo_Popurri_RandomForest') %>%
show_best(metric = 'sensitivity', n=50) %>%
filter(row_number()==1) %>%
select(mtry,trees,min_n,.config )
Modelo_Popurri_RandomForest_final <- Wokflows_set_map %>%
extract_workflow('Modelo_Popurri_RandomForest') %>%
finalize_workflow(Modelo_Popurri_RandomForest_hyperparametros) %>%
last_fit(Split_datos_strat, metrics = Modelo_Metricas )Una vez acabado todo, podemos ver todo lo necesario usando unos pocos verbos como collect y las lógicas del tidyverse.
map2(
list(
Modelo_Escalas_RandomForest_final,
Modelo_Cognicion_RandomForest_final,
Modelo_Bioquimica_RandomForest_final,
Modelo_Popurri_RandomForest_final),
list('Modelo_Escala','Modelo_Cognicion','Modelo_Bioquimica','Modelo_Popurri'),
~ ..1 %>%
collect_metrics() %>%
mutate('Modelo'= ..2) %>%
select(.metric,.estimate,Modelo )) %>%
bind_rows() %>%
pivot_wider(names_from = Modelo, values_from = .estimate) %>%
knitr::kable()| .metric | Modelo_Escala | Modelo_Cognicion | Modelo_Bioquimica | Modelo_Popurri |
|---|---|---|---|---|
| accuracy | 0.6666667 | 0.6091954 | 0.7241379 | 0.6436782 |
| j_index | 0.2922705 | -0.0265700 | -0.0048309 | 0.1400966 |
| precision | 0.8703704 | 0.7868852 | 0.7922078 | 0.8275862 |
| sensitivity | 0.6811594 | 0.6956522 | 0.8840580 | 0.6956522 |
| specificity | 0.6111111 | 0.2777778 | 0.1111111 | 0.4444444 |
| f_meas | 0.7642276 | 0.7384615 | 0.8356164 | 0.7559055 |
| recall | 0.6811594 | 0.6956522 | 0.8840580 | 0.6956522 |
| mcc | 0.2440012 | -0.0235126 | -0.0061354 | 0.1203859 |
| roc_auc | 0.7302738 | 0.4927536 | 0.4919485 | 0.6344605 |
ROC_CURVE_training_todos_los_posibles_modelos <- Wokflows_set_map %>%
collect_predictions() %>%
group_split(wflow_id) %>%
set_names(list('Escala','Cognicion','Bioquimica','Popurri')) %>%
map(~ .x %>% roc_curve(truth=Var_respuesta , .pred_0)) %>%
map2(., list('Escala','Cognicion','Bioquimica','Popurri'),
~ .x %>% mutate(model= .y)) %>%
bind_rows()
Plot_ROC_CURVE_training_todos_los_posibles_modelos <- ROC_CURVE_training_todos_los_posibles_modelos %>%
ggplot(aes( x = 1- specificity ,y= sensitivity, color= model ) ) +
geom_path() +
geom_abline(intercept=0, slope=1, linetype=3)
ROC_CURVE_training_modelos_finales <- list(
list('Esc' ,Modelo_Escalas_RandomForest_hyperparametros$.config),
list('Cog' ,Modelo_Cognicion_RandomForest_hyperparametros$.config),
list('Bio' ,Modelo_Bioquimica_RandomForest_hyperparametros$.config),
list('Pop' ,Modelo_Popurri_RandomForest_hyperparametros$.config)) %>%
map(~ Wokflows_set_map %>%
collect_predictions() %>%
filter(
str_detect(wflow_id, .x[[1]]) &
.config== .x[[2]]
)) %>%
map(~ .x %>%
mutate(Model= str_extract(wflow_id, '^.{2}') ) %>%
select(Model,Var_respuesta, .pred_0, .pred_1, .pred_class)
) %>%
set_names(list('Escala','Cognicion','Bioquimica','Popurri') ) %>%
map(~ .x %>% roc_curve(truth = Var_respuesta, .pred_0)) %>%
map2(., list('Escala','Cognicion','Bioquimica','Popurri'),
~ .x %>% mutate(model= .y)) %>%
bind_rows()
Plot_ROC_CURVE_training_modelos_finales <- ROC_CURVE_training_modelos_finales %>%
ggplot(aes( x = 1- specificity ,y= sensitivity, color= model ) ) +
geom_path() +
geom_abline(intercept=0, slope=1, linetype=3)
ROC_CURVE_testing_modelos_finales <-list(
Modelo_Escalas_RandomForest_final,
Modelo_Cognicion_RandomForest_final,
Modelo_Bioquimica_RandomForest_final,
Modelo_Popurri_RandomForest_final
) %>%
map(~.x %>%
collect_predictions() %>%
roc_curve(Var_respuesta, .pred_0)) %>%
map2(., list('Escala','Cognicion','Bioquimica','Popurri'),
~ .x %>% mutate(model= .y)) %>%
bind_rows()
Plot_ROC_CURVE_testing_modelos_finales <- ROC_CURVE_testing_modelos_finales %>%
ggplot(aes(x = 1- specificity ,y = sensitivity, color = model ) ) +
geom_path()+
geom_abline(intercept=0, slope=1, linetype=3)ggpubr::ggarrange( Plot_ROC_CURVE_training_todos_los_posibles_modelos,
Plot_ROC_CURVE_training_modelos_finales,
Plot_ROC_CURVE_testing_modelos_finales ,
nrow=1 ,
common.legend = T,
legend="bottom",
labels = c('Todos_training','Training','Test')) Los valores de Shapley (shap values), que deben su nombre al premio Nobel Lloyd Shapley, son un concepto de la teoría de juegos cooperativos que ha acabado encontrado aplicaciones en el aprendizaje automático.
En teoría de juegos, los valores de Shapley son la distribución equitativa de la contribución de cada jugador en un juego cooperativo entre todos los jugadores. En el contexto del aprendizaje automático, las características o variables se consideran jugadores, y el juego consiste en cuánto contribuye cada característica a la predicción de un modelo.
Es decir, los valores de Shapley son sencillamente las marginales de cada variable dentro del modelo. Compilarlos es computacionalmente extenso, ya que hay que hacer cada combinación de valores para una fila y evaluar de cuanto es su marginal. El valor final llamado Shap value es promedio de todas las aportaciones de una variable marginal.
Aparte de su extenso coste computacional, todo shap value asume que está utilizando un modelo de caja negra. de manera que si el modelo no está bien especificado en tema de interacciones, etc. es posible que estén dando indicaciones erróneas.
Actualmente se pueden compilar en tidymodels utilizando la libreria agua, al ser una extensión de Tidymodels del entorno de h2o. El problema es que esto hace que el motor de los modelos deba ser a la fuerza h2o.
Una propuesta sería buscar un modelo bien calibrado con Tidymodels y pasarlo por h2o para hacer un diagnostico más eficaz.
h2o_start()
modelo_h20 <- rand_forest() %>%
set_engine("h2o", max_runtime_secs = 20) %>%
set_mode('classification') %>%
fit(Var_respuesta ~ ., data = Receta_modelo_Bioquimica %>% bake(new_data = NULL) )
entrada_h20 <- h2o::as.h2o(
Receta_modelo_Bioquimica %>% bake(new_data = NULL) %>%
dplyr::rename(.outcome = Var_respuesta))##
|
| | 0%
|
|======================================================================| 100%
Para cada predicción, los valores de Shapley representan la contribución de cada característica a la diferencia entre la predicción del modelo y la predicción media. En otras palabras, indican cuánto añade o resta cada característica a la predicción media.
Si una predicción media es de 0.5, y el shap value es de 0.1, indicará que esa variable sola podría llegar a tener un impacto tal que haga que la predicción llege a 0.6.